{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Logistic regression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.导入需要的库/模块"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "from logistic_regression import LogisticRegression"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.加载MNIST数据集\n",
    "\n",
    "在** ./data/code/transform.c **中，我们对数据集进行了转换和拆分为训练集(training)、验证集(validation)和测试集(testing)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Loading training data\n",
      "Loading testing data\n",
      "Loading validating data\n"
     ]
    }
   ],
   "source": [
    "print(\"Loading training data\")\n",
    "pd_train = pd.read_csv(\"../data/training\", header=None)\n",
    "\n",
    "print(\"Loading testing data\")\n",
    "pd_test = pd.read_csv(\"../data/testing\", header=None)\n",
    "\n",
    "print(\"Loading validating data\")\n",
    "pd_validate = pd.read_csv(\"../data/validation\", header=None)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "查看数据集的维度信息："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(50000, 785)\n",
      "(10000, 785)\n",
      "(10000, 785)\n"
     ]
    }
   ],
   "source": [
    "print(pd_train.shape)\n",
    "print(pd_test.shape)\n",
    "print(pd_validate.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们可以看到训练集有50000个样本，验证集和测试集有10000个样本，每一个样本都是1+28*28=785维的特征向量"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.训练模型\n",
    "\n",
    "MNIST一共类10类(数字0-9)，因此我们需要将Logistic regression推广到多分类问题上。\n",
    "\n",
    "基本思路：将多分类任务拆若干个二分类任务。学习时分别训练多个分类器，测试时对这个分类器的预测结果进行集成以获得最终的分类结果。\n",
    "\n",
    "经典的拆分策略有：\n",
    "\n",
    "- 一对一(OvO)\n",
    "- 一对多(OvR)\n",
    "- 多对多(MvM)\n",
    "\n",
    "我们将使用一对多的拆分方法，每次将一个类的样例作为正例、所有其他类的样例作为反例来训练10个分类器。在测试时选择sigmoid函数值最大的，所对应的类标号作为最终的分类结果。\n",
    "\n",
    "**baseline:**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "------------- Training model for number: 0 -------------\n",
      "\n",
      "step 0: 0.702554\n",
      "step 10: 0.098281\n",
      "step 20: 0.079234\n",
      "step 30: 0.069833\n",
      "step 40: 0.063964\n",
      "step 50: 0.059857\n",
      "step 60: 0.056778\n",
      "step 70: 0.054360\n",
      "step 80: 0.052397\n",
      "step 90: 0.050761\n",
      "Final training error: 0.049501\n",
      "\n",
      "------------- Training model for number: 1 -------------\n",
      "\n",
      "step 0: 0.519057\n",
      "step 10: 0.083981\n",
      "step 20: 0.068974\n",
      "step 30: 0.061465\n",
      "step 40: 0.056826\n",
      "step 50: 0.053625\n",
      "step 60: 0.051261\n",
      "step 70: 0.049429\n",
      "step 80: 0.047960\n",
      "step 90: 0.046750\n",
      "Final training error: 0.045827\n",
      "\n",
      "------------- Training model for number: 2 -------------\n",
      "\n",
      "step 0: 0.718578\n",
      "step 10: 0.165012\n",
      "step 20: 0.137207\n",
      "step 30: 0.123960\n",
      "step 40: 0.115907\n",
      "step 50: 0.110379\n",
      "step 60: 0.106295\n",
      "step 70: 0.103122\n",
      "step 80: 0.100566\n",
      "step 90: 0.098451\n",
      "Final training error: 0.096830\n",
      "\n",
      "------------- Training model for number: 3 -------------\n",
      "\n",
      "step 0: 0.727284\n",
      "step 10: 0.173113\n",
      "step 20: 0.148019\n",
      "step 30: 0.135902\n",
      "step 40: 0.128520\n",
      "step 50: 0.123455\n",
      "step 60: 0.119714\n",
      "step 70: 0.116806\n",
      "step 80: 0.114460\n",
      "step 90: 0.112512\n",
      "Final training error: 0.111012\n",
      "\n",
      "------------- Training model for number: 4 -------------\n",
      "\n",
      "step 0: 0.639896\n",
      "step 10: 0.150056\n",
      "step 20: 0.123644\n",
      "step 30: 0.110424\n",
      "step 40: 0.102179\n",
      "step 50: 0.096428\n",
      "step 60: 0.092137\n",
      "step 70: 0.088785\n",
      "step 80: 0.086079\n",
      "step 90: 0.083840\n",
      "Final training error: 0.082126\n",
      "\n",
      "------------- Training model for number: 5 -------------\n",
      "\n",
      "step 0: 0.635689\n",
      "step 10: 0.198403\n",
      "step 20: 0.170377\n",
      "step 30: 0.155645\n",
      "step 40: 0.146322\n",
      "step 50: 0.139796\n",
      "step 60: 0.134926\n",
      "step 70: 0.131124\n",
      "step 80: 0.128054\n",
      "step 90: 0.125513\n",
      "Final training error: 0.123564\n",
      "\n",
      "------------- Training model for number: 6 -------------\n",
      "\n",
      "step 0: 0.686778\n",
      "step 10: 0.122983\n",
      "step 20: 0.096807\n",
      "step 30: 0.084995\n",
      "step 40: 0.078023\n",
      "step 50: 0.073325\n",
      "step 60: 0.069898\n",
      "step 70: 0.067262\n",
      "step 80: 0.065158\n",
      "step 90: 0.063429\n",
      "Final training error: 0.062113\n",
      "\n",
      "------------- Training model for number: 7 -------------\n",
      "\n",
      "step 0: 0.631441\n",
      "step 10: 0.119841\n",
      "step 20: 0.098167\n",
      "step 30: 0.088150\n",
      "step 40: 0.082173\n",
      "step 50: 0.078115\n",
      "step 60: 0.075136\n",
      "step 70: 0.072830\n",
      "step 80: 0.070978\n",
      "step 90: 0.069446\n",
      "Final training error: 0.068273\n",
      "\n",
      "------------- Training model for number: 8 -------------\n",
      "\n",
      "step 0: 0.785142\n",
      "step 10: 0.249144\n",
      "step 20: 0.218743\n",
      "step 30: 0.203384\n",
      "step 40: 0.193617\n",
      "step 50: 0.186632\n",
      "step 60: 0.181265\n",
      "step 70: 0.176934\n",
      "step 80: 0.173315\n",
      "step 90: 0.170211\n",
      "Final training error: 0.167752\n",
      "\n",
      "------------- Training model for number: 9 -------------\n",
      "\n",
      "step 0: 0.697977\n",
      "step 10: 0.209030\n",
      "step 20: 0.183557\n",
      "step 30: 0.169284\n",
      "step 40: 0.159904\n",
      "step 50: 0.153187\n",
      "step 60: 0.148099\n",
      "step 70: 0.144088\n",
      "step 80: 0.140828\n",
      "step 90: 0.138115\n",
      "Final training error: 0.136028\n",
      "\n",
      "Evaluate model on validating data\n",
      "\n",
      "On validating data: 8966/10000\n",
      "\n"
     ]
    }
   ],
   "source": [
    "PREDICT_NUMBERS = 10\n",
    "\n",
    "train_X = pd_train.drop([0], axis=1).values / 255.0\n",
    "predictors = list()\n",
    "for number in range(0, PREDICT_NUMBERS):\n",
    "    train_y = np.array(list(map(int, pd_train[0].values == number)))\n",
    "    predictors.append(LogisticRegression())\n",
    "    print(\"------------- Training model for number: %d -------------\\n\" % number)\n",
    "    predictors[number].fit(train_X, train_y, 0.5, batch_num=50000, verbose=True)\n",
    "\n",
    "print(\"Evaluate model on validating data\")\n",
    "validate_X = pd_validate.drop([0], axis=1).values / 255.\n",
    "tot_right_num = 0\n",
    "for tx, ty in zip(validate_X, pd_validate[0].values):\n",
    "    pred_y = list()\n",
    "    for number in range(0, PREDICT_NUMBERS):\n",
    "        pred_y.append(predictors[number].predict(np.reshape(tx, (1, -1))))\n",
    "\n",
    "    if np.argmax(pred_y) == ty:\n",
    "        tot_right_num += 1\n",
    "\n",
    "print(\"\\nOn validating data: %s/%s\\n\" % (tot_right_num, validate_X.shape[0]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们看到该模型在验证集上的准确率有89.6%，这是个不错的结果。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4.参数调优\n",
    "\n",
    "#### 学习率\n",
    "\n",
    "我们设置7个不同的学习率，各训练一个模型，然后选出在验证集上表现最好的作为最终的学习率："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "------------- Training model with alpha 0.001 -------------\n",
      "On validating data: 9175/10000\n",
      "On training data: 45497/50000\n",
      "------------- Training model with alpha 0.003 -------------\n",
      "On validating data: 9212/10000\n",
      "On training data: 45878/50000\n",
      "------------- Training model with alpha 0.010 -------------\n",
      "On validating data: 9237/10000\n",
      "On training data: 46144/50000\n",
      "------------- Training model with alpha 0.030 -------------\n",
      "On validating data: 9231/10000\n",
      "On training data: 46245/50000\n",
      "------------- Training model with alpha 0.100 -------------\n",
      "On validating data: 9219/10000\n",
      "On training data: 46120/50000\n",
      "------------- Training model with alpha 0.300 -------------\n",
      "On validating data: 9116/10000\n",
      "On training data: 45576/50000\n",
      "------------- Training model with alpha 1.000 -------------\n",
      "On validating data: 8770/10000\n",
      "On training data: 43720/50000\n"
     ]
    }
   ],
   "source": [
    "accuracy_tra = []\n",
    "accuracy_val = []\n",
    "\n",
    "for alpha in [0.001, 0.003, 0.01, 0.03, 0.1, 0.3, 1]:\n",
    "    predictors = list()\n",
    "    print(\"------------- Training model with alpha %.3f -------------\" % alpha)\n",
    "    for number in range(0, PREDICT_NUMBERS):\n",
    "        train_y = np.array(list(map(int, pd_train[0].values == number)))\n",
    "        predictors.append(LogisticRegression())\n",
    "    #     print(\"------------- Training model for number: %d -------------\\n\" % number)\n",
    "        predictors[number].fit(train_X, train_y, alpha, batch_num=10, verbose=False)\n",
    "\n",
    "    # print(\"Evaluate model on validating data\")\n",
    "    tot_right_num_val = 0\n",
    "    for tx, ty in zip(validate_X, pd_validate[0].values):\n",
    "        pred_y = list()\n",
    "        for number in range(0, PREDICT_NUMBERS):\n",
    "            pred_y.append(predictors[number].predict(np.reshape(tx, (1, -1))))\n",
    "\n",
    "        if np.argmax(pred_y) == ty:\n",
    "            tot_right_num_val += 1\n",
    "\n",
    "    print(\"On validating data: %s/%s\" % (tot_right_num_val, validate_X.shape[0]))\n",
    "    \n",
    "    tot_right_num_tra = 0\n",
    "    for tx, ty in zip(train_X, pd_train[0].values):\n",
    "        pred_y = list()\n",
    "        for number in range(0, PREDICT_NUMBERS):\n",
    "            pred_y.append(predictors[number].predict(np.reshape(tx, (1, -1))))\n",
    "\n",
    "        if np.argmax(pred_y) == ty:\n",
    "            tot_right_num_tra += 1\n",
    "\n",
    "    print(\"On training data: %s/%s\" % (tot_right_num_tra, train_X.shape[0]))\n",
    "    \n",
    "    accuracy_tra.append(tot_right_num_tra / 50000)\n",
    "    accuracy_val.append(tot_right_num_val / 10000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x90c8240>]"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX8AAAEACAYAAABbMHZzAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmcjXX/x/HXZxZrZS8iJFlG9kJpGfsod6vuUrmjQlkj\nYysapSyDEC2qOy1Kt1ZtvxJGpUQZDGYsLYqSqOzLLN/fH9epJmHOMHPOmTnv5+NxHs51ru/3nM+5\nHuNzXed7fRdzziEiIuElItgBiIhI4Cn5i4iEISV/EZEwpOQvIhKGlPxFRMKQkr+ISBjyK/mbWZyZ\npZnZBjMbepT9pc3sdTNbZWZLzSzG93pRM/vCzJLNLMXM7s/rLyAiIrlnOfXzN7MIYAPQBvgRWA7c\n6JxLy1ZmArDHOfegmdUGZjjn2vr2lXDO7TezSGAJ0N85tyx/vo6IiPjDnyv/ZsBG59xm51w6MAe4\n6ogyMcBCAOfceqC6mVXwbe/3lSkKRAEaVSYiEmT+JP/KwA/Ztrf4XstuFXAtgJk1A6oCVXzbEWaW\nDGwD5jvnlp9s0CIicnLy6obvOKCMma0A+gDJQCaAcy7LOdcY72TQ/I/7ASIiEjxRfpTZincl/4cq\nvtf+5JzbA9z2x7aZfQt8c0SZ3Wa2CIgD1h35IWam5iARkVxyztmJ1PPnyn85UNPMqplZEeBGYF72\nAmZWysyifc97AIudc3vNrLyZlfK9XhxoB6RxDM45PZzj/vvvD3oMofDQcdCx0LE4/uNk5Hjl75zL\nNLO+wId4J4tnnHOpZtbL2+1mAnWB58wsC1gL3O6rXsn3eoSv7ivOufdOKmIRETlp/jT74Jz7P6D2\nEa89me350iP3+15PAZqcZIwiIpLHNMI3BMXGxgY7hJCg4/AXHYu/6FjkjRwHeQWKmblQiUVEpCAw\nM1w+3vAVEZFCRslfRCQMKfmLiIQhJX8RkTCk5C8iEoaU/EVEwpCSv4hIGFLyFxEJQ0r+IiJhSMlf\nRCQMKfmLiIQhJX8RkTBUYJL/li3QvDl07w6ffw6aA05E5MQViOT/++/QsaP3iImB//wHGjSAadPg\nt9+CHZ2ISMET8lM6HzoEcXFesp8yBcy8q/7Fi2HmTHjvPfjXv6BnT7j4Ym+/iEg4OJkpnUM6+Wdl\nwU03QWYmzJkDkZH/rLdjB7zwgnciAOjRw/tlUL58AIIWEQmiQjWf/4qfVvD2+rc5lHGI+Hj48Ucv\nuR8t8YOX5AcOhHXr4KmnYOVKqFkTunSBRYt0b0BE5GhC7sq/4+yOfPvbt3y/82eKfNeJx/v+m2sb\ntqdoVFG/3+u33+DFF71fA4cOwR13QLducPrp+Re/iEigFZor//TMdJZ8v4RBpT+l9Ox13HNjc55I\nmUjFSRXp+kZX3kp7ix37d+T4PmXKQL9+sHo1PP88pKVBrVpw/fUwf77XnCQiEs5C6sr/8x8+55Y5\nvdg9YRULFkD9+t6+n/b8xGupr/FG2ht8+eOXnFb0NBpXbEyjio1oXLExjSs1plqpathx7vbu2gUv\nveT9Gti1y/s10L07VKoUoC8oIpLHCs0N37GfjOW5136i37lT6d376OWcc3z7+7ck/5RM8rZkVm5b\nSfK2ZA6kH6BRxUZ/OyHUKV+HqIioI+rDV195J4G5cyE21usp1L79se8riIiEokKT/Nu/0J6tb97F\ntN5X07p17upv37f9HyeELbu3UK9Cvb+dEBqc0YAS0SUA2LPH60U0cyZs3+7dOO7TB6Kj8+ELiojk\nsUKT/E99+FRKz/qOhe+WpWbNk3/PPYf2sPrn1X87IaT+kkr10tX/dkJoXLExm9PKMWwYbNsGM2bA\nJZec/OeLiOSnQpP8Gz3eiHUDktm1C4oVy5/POZx5mNRfUv88GfxxYihVtBSdav2LpnsSGD2kAq1a\nwYQJcMYZ+ROHiMjJOpnkH5VzkcA5v0Is28rmX+IHKBJZhIYVG9KwYkNu5VYAslwW3/z2DdOXTWfY\nDzEMfH4YO9/rx3nnFeH+++Guu3Q/QEQKl5Dq6lkzshVVqwb+cyMsgpplazIlbgqfdP+EJT8uZN5Z\n9bj/pXm8+prjggtg6dLAxyUikl/8Sv5mFmdmaWa2wcyGHmV/aTN73cxWmdlSM4vxvV7FzBaa2Voz\nSzGz/sf7nFJ7mwcl+WdXp3wd3r3pXR7t+CiPbRhGdPf2/LtfCtdd53UP3ZHzMAMRkZCXY/I3swhg\nOtABqAd0MbM6RxQbASQ75xoCtwLTfK9nAIOcc/WAC4E+R6n7p1+3lqFatdx/ifwQVzOOVXeu4sra\nVzJ5Rxs6Tu9N5Gm/UK+e1ztIA8VEpCDz58q/GbDRObfZOZcOzAGuOqJMDLAQwDm3HqhuZhWcc9uc\ncyt9r+8FUoHKx/qgrd8XCfqVf3bRkdH0a96PtL5plCwWzeuVYuj62CPMeuEwLVp44wVERAoif5J/\nZeCHbNtb+GcCXwVcC2BmzYCqQJXsBcysOtAI+OJYH/T994RU8v9D2eJlmdpxKh93+5i1Bz9k57/r\nc1G3d7j8Ckfv3lpTQEQKnry64TsOKGNmK4A+QDKQ+cdOMzsFeBUY4PsFcFTff0/INPscTd0KdXn/\n5veZEjeFD1w8MQ934LeotdStC88+q6YgESk4cuznb2YtgATnXJxvexjgnHPjj1PnW6C+c26vmUUB\n7wDvO+emHqeOK1r0fvr3hxIlIDY2ltjY2BP4SoGRnpnO418+zpiPxxBb4Xo2PT2aEpRnxgxo2DDY\n0YlIYZSUlERSUtKf26NHj86/QV5mFgmsB9oAPwHLgC7OudRsZUoB+51z6WbWA2jpnOvm2/c8sMM5\nNyiHz3ElSjj27i1Yq3Ht3L+ThKQE5qydQ+uoe1k0oQ833RDN6NFQqlSwoxORwixfp3R2zmUCfYEP\ngbXAHOdcqpn1MrOevmJ1gTVmlorXK2iAL7CWwM1AazNLNrMVZhZ3rM+qVq1gJX6AciXK8ejlj7K4\n22J2Vfg/Sg2rT2rmu9SNccyercVkRCQ0hdT0Dh06OP7v/4IdyYlzzvH+pvcZ9MEgSlOdXXMnUzEi\nhhkzvIXnRUTyUqFZzCUUe/rkhplx+bmXk3JXCl3O78iOTpeRFdePSzvsZMgQ2HvMW90iIoEVUsk/\nlHv65EZ0ZDQDWgwgtU8q9es76FuXxQenUScmnblz1RQkIsEXUsm/oF/5H6l8ifJMv3w6i7svotQF\n7xDZrwFDZr5Phw6wfn2woxORcKbkHwD1Tq/HB7d8wIyrEom6cgA/XNqR5p1Sufde2L8/2NGJSDhS\n8g8QM6NTrU6s7b2Gnm3aE3nHpby6rz+1G/3KW2+pKUhEAiukevv89pujdOlgRxIYO/bvYNSiUby8\n6lWKLB1Jk6w7mTEtmho1gh2ZiBQUhWYlr337HCVKBDuSwFqzfQ0D3h/Ims1bOfDmZO65Ko6hQ/N3\nQRsRKRwKTfJPT3dEhdTaYoHhnOOdDe8w4L17OLDlXKIXTeKJMXW4/PJgRyYioazQJP+sLFfgRvjm\npcOZh5m+bDoPLBwLKTfTMmMUj00qW2i6wIpI3io0g7zCOfGDt77woAsHsfHudfy7yyEW169LvW4z\nGPNwBocPBzs6ESlMQurKP1RiCRWrf17NnW8OJOXbbZRe+gjP3teetm2DHZWIhIpC0+wTKrGEEucc\n89bP484372HXprpcdmgiT4+rTeVjrocmIuFCyT8MHMo4xOQlj/LgonG4VV0Z0XIUwwaUITo62JGJ\nSLAUmjZ/ObaiUUUZftlgvotfx1Wd9zPmtzpU6/wYCxZlBDs0ESmAdOVfQK38aRVdZw9k/ZbtXLr/\nEV58oB0VKwY7KhEJJDX7hCnnHK+seou73hjM/s0xxDecSEK/WmE5VkIkHCn5h7lDGYe49+1pTPtq\nPKW/u5UXeoykw2VhMk+GSBhT8hcAtu35mZufGUnST/NofiCB1+69g0pn6GeASGGl5C9/8+mmlXR5\nbiA/7dpB33OmMKlvGyIjgx2ViOQ1JX/5B+ccUz54gxFJgym2qz5P3zCR62LPDXZYIpKHlPzlmA6k\nH6TbE1OZuzWR+pndeWvQfVSvVCrYYYlIHlA/fzmm4tHFeKXfUFL7rSEj6jfOmVKb/0x9kvSMzGCH\nJiJBpCv/MPPSomTufPNuMqN/55EOj9CzXetghyQiJ0jNPpIrmZmOfo+9zszv4qlWrCFzeyTSpHrN\nYIclIrmk5C8nZOvPB7l2whSWR0+kY8XbmN3rXkoX1/0AkYJCbf5yQiqfUYwvJg3jnctT+CJlJ2c8\nWIeEt58iM0v3A0QKO135CwCZmXDfjK+YtHYgpc/YzX9vmEynerofIBLK1OwjeWb7dscND7zKJ0WH\ncV7FOrz4nwmcd0a9YIclIkeh5C957uPPDnHLlCfYVvNhrjj3X8zo/ABnnnpmsMMSkWzyvc3fzOLM\nLM3MNpjZ0KPsL21mr5vZKjNbamYx2fY9Y2Y/m9nqEwlQguPSi4ry3ZwBTKy+ng/fKkeNifUZ/N4o\n9hzaE+zQRCQP5Jj8zSwCmA50AOoBXcyszhHFRgDJzrmGwK3AtGz7nvXVlQImIgL69yzNlmfH02VX\nMtNf+I4qE2oxY9njpGemBzs8ETkJ/lz5NwM2Ouc2O+fSgTnAVUeUiQEWAjjn1gPVzayCb/tT4Le8\nC1kCrUwZeHZKVZaPeJ5zl73HsOdeo+bk+ryZ9iZqqhMpmPxJ/pWBH7Jtb/G9lt0q4FoAM2sGVAWq\n5EWAEjrq14fl8xrz9KXz2f/6FLo9P4oWT17K0i1Lgx2aiORSXk32Pg6YamYrgBQgGch1Z/GEhIQ/\nn8fGxhIbG5tH4UleMYMbbjA6dYrjobHtmDb7eeJ+vp42tVswvv1YapbVSGGR/JKUlERSUlKevFeO\nvX3MrAWQ4JyL820PA5xzbvxx6nwL1HfO7fVtVwPeds41OE4d9fYpgL7+GvoN2s+yiCmknz+Zbk1u\nZuRlIylfonywQxMp9PK7t89yoKaZVTOzIsCNwLwjAihlZtG+5z2AxX8k/j+K+B5SyJxzDrz3Vgle\n6DmC8q+kMu/tLGpNq8O4T8dxIP1AsMMTkWPIMfk75zKBvsCHwFpgjnMu1cx6mVlPX7G6wBozS8Xr\n2TPgj/pm9hLwGVDLzL43s+55/SUk+Dp2hNQvK3BX9UdxT3/GrA+/pNajtXlu5XOaLkIkBGmQl+S5\nH3+EoUPhg3WfUfrfgylRaj+J7RJpd067YIcmUqhohK+EpCVLoE9fR/o5b7C3xTDqVDybCW0n0LBi\nw2CHJlIoaFZPCUktW8JXXxr92l7LgUlrObjqStq/0IFub3Zjy+4twQ5PJKwp+Uu+ioyEO++E1DXR\n1NvXh6ypG9j+dRUaPtGQEQtGsOvgrmCHKBKWlPwlIMqVg8ceg/nvnMaeN8dQed4qUr79mVrTazHt\ni2kczjwc7BBFwora/CXgnIOXX4YhQ6BJXAp7LxrCD/s2MbbNWK6rex1m6hUs4g/d8JUCac8eeOgh\nePppuDb+I5aViqd4dDEmtptIy6otgx2eSMhT8pcCbcMGuPtu+PqbLK4aOZtXfrmPppWaMrbNWGqX\nrx3s8ERClnr7SIFWqxa8+y5MTIzgtVFdafjxemqVaMHFz15Mn3f7sH3f9mCHKFLoKPlLSDCDf/0L\n1q6FFucX4+nuQ7h1TxrmihAzI4YxH49hf/r+YIcpUmgo+UtIKVYMRoyA5GTYsqEcb/d5hIQqX7Bm\n+xpqPVqLZ1Y8o+kiRPKA2vwlpCUlQf/+UKEC3HH/Fzz+dTy/HviVCe0m0LFmR/UMkrCmG75SqGVk\nwOOPwwMPwM23OJp1fZsHPx/KmaeeSWK7RJpUahLsEEWCQjd8pVCLioJ+/WDdOti/z7jniisZVCyF\nznX/TaeXOnHL67ew+ffNwQ5TpEDRlb8UOF9+6Z0MsrJg3CN7WHx4Eo8ue5TbGt3GiEtGUKZ4mWCH\nKBIQuvKXsHL++d6Mob17w82dT+WHFxJIun4Nuw7tovb02jzy+SMcyjgU7DBFQpqSvxRIERFw662Q\nlgZlykDrZpWo981M5t+cxMLvFlJ3Rl3mrJlDlssKdqgiIUnNPlIopKbCgAHeQjLTpkFEjSQGfziY\nCIsgsV0il1W/LNghiuQ59fYRwZsw7s03YdAgr2kocWIWn+9+hRELR1D/9PqMazuOmAoxwQ5TJM+o\nzV8Eb5TwNdd4vYLq14fzm0aw6c0urLw9jdjqscTOiqXn2z35ac9PwQ5VJOiU/KXQKV4cRo2Cr76C\nVaugcYOi1Ph5EGl91lOqaCnOe/w8EpIS2Ht4b7BDFQkaNftIoffRR94o4bPOgqlToVjF77hv4X0s\n/HYh9192P7c3uZ2oiKhghymSa2r2ETmOtm29XwBxcXDJJTDjoeo81vZF3u7yNq+sfYX6j9dn3vp5\n6OJDwomSv4SF6GgYOBDWrIGdO6FuXVj7UVPm37KASe0nMWLBCGKfi2X51uXBDlUkINTsI2Hpiy+8\nUcJRUTB9OjRolMFzK59jVNIoLql6CQ+3eZgaZWoEO0yR41Kzj0guNW8OS5fC7bfD5ZdDn7uiuKrq\n7Wzou4F6FepxwVMXMPD/BrJz/85ghyqSL5T8JWxFRHjJPy3NW0cgJgZmPVWS4S1Hsq73Og5lHqLO\njDpMWDKBgxkHgx2uSJ5Ss4+Iz5o1Xq+gnTvh0Ufh0kth/Y71DFswjBU/rWBMqzHc3OBmIkzXTBIa\nNMJXJI84B6++CoMHw0UXQWIiVKkCn37/KfHz4zmUcYjEdom0qdEm2KGK5H+bv5nFmVmamW0ws6FH\n2V/azF43s1VmttTMYvytKxJKzOD6671RwjVrQqNGMHYsXHDGxXx222eMuGQEvd7pRcfZHUn5OSXY\n4YqcsByTv5lFANOBDkA9oIuZ1Tmi2Agg2TnXELgVmJaLuiIhp2RJePBBWLbM6xlUrx68+67ROaYz\n6/qso2PNjrR9oS23v3U7W3dvDXa4Irnmz5V/M2Cjc26zcy4dmANcdUSZGGAhgHNuPVDdzCr4WVck\nZNWo4U0WN3063HMPXHEFbP6mCP2b92dD3w2cXvJ0GjzRgPsW3sfuQ7uDHa6I3/xJ/pWBH7Jtb/G9\nlt0q4FoAM2sGVAWq+FlXJOTFxUFKCsTGwoUXwvDhEJlRirFtx7Ky10q27N5CrUdrMWPZDNIz04Md\nrkiO8mpCk3HAVDNbAaQAyUBmbt8kISHhz+exsbHExsbmUXgiJ69IEYiPh5tvhqFDvVHCEybAjTee\nxayrZ7Fy20qGzB/C1C+mMq7tOK6pcw1mJ3QvTuSokpKSSEpKypP3yrG3j5m1ABKcc3G+7WGAc86N\nP06db4H6wHn+1lVvHylolizxRgmfcorXNbRhQ+/1D7/+kPj58ZxS5BQS2yVy0VkXBTdQKbTyu7fP\ncqCmmVUzsyLAjcC8IwIoZWbRvuc9gMXOub3+1BUpqFq2hOXLvV8C7dtDnz7w66/Q/pz2rOi5gp5N\nenLjqzfS+X+d2bhzY7DDFfmbHJO/cy4T6At8CKwF5jjnUs2sl5n19BWrC6wxs1S8nj0Djlc377+G\nSHBERkKvXt4ykuA1BT35JOAiubXRrazvu57zzzyfC5+5kH7v9eOXfb8ENV6RP2iQl0geWrXKawra\nu9drCmrZ0nv9l32/MObjMcxOmc2gCwdxd4u7KRFdIrjBSoGnEb4iIcQ5mDPHuzncqpV3U7hSJW/f\npl83MXzBcJZuWcqDrR6ka4OuREZEBjdgKbA0q6dICDGDLl28CeOqVPHWE05MhMOHoWbZmsy9fi5z\nr5/L0yuepsnMJnyw6YNghyxhSFf+Ivls40a4+274+mtvGckOHbzXnXO8mfYmwxYMo1qpakxoN4FG\nFRsFN1gpUNTsI1IAvPOOdxI47zyYPNkbPQyQnpnOUyue4oHFD9D+nPaMaT2GqqWqBjdYKRDU7CNS\nAHTqBGvXegvJNGsGI0fC/v0QHRlN7wt6s6HfBqqVqkbjJxszdP5Qfj/4e7BDlkJMyV8kgIoW9aaG\nWLkSNm3yuobOnevdJD6t6Gk82PpBUu5KYeeBndSeXpupS6dyOPNwsMOWQkjNPiJBtHix1zW0fHmv\na2i9en/tW7N9DUM/GkrajjTGthnL9THXa7oI+Ru1+YsUYBkZ3sCw0aPhppsgIQFKl/5r/8JvFxI/\nP56oiCgmtpvIJdUuCVqsElrU5i9SgEVFeVNDrFsHBw5AnTrwzDOQleXtb312a5b3WE7/Zv3p+kZX\nrp5zNWk70oIbtBR4uvIXCTFffeU1BWVkeE1BzZv/te9gxkGmL5vO+CXj6Vy3MwmxCZxxyhnBC1aC\nSlf+IoVI06bw6afQty9ccw3cdhv8/LO3r1hUMQZfNJi0PmkUjy5OzGMxPLD4AfYd3hfcoKXAUfIX\nCUEREfCf/3ijhMuV88YGPPIIpPvWiSlXohyTO0zmyx5fkrojlVrTa/HUV0+RkZUR3MClwFCzj0gB\nkJYGAwbAli0wbRq0afP3/cu3Lid+fjy/7P+F8W3Hc8W5V6hnUBhQbx+RMOAcvPUWDBoETZrApElQ\nrVr2/Y53N77LkPlDOOOUM0hsl8j5Z54fvIAl36nNXyQMmMHVV3ujhBs08E4Ao0d7PYS8/UanWp1Y\nfddqbjrvJq58+Upueu0mvvv9u6DGLaFJyV+kgCleHEaNghUrYM0aiImBN97wfhkAREVE0aNpDzb0\n20DtcrVpOrMpgz8czK8Hfg1u4BJS1OwjUsAtWAD9+0Plyt79gDp1/r5/295tJCQl8Hrq6wxtOZS+\nzfpSNKpocIKVPKVmH5Ew1qaNN1fQ5ZfDJZfA4MGwe/df+yueUpEnOj3B4m6L+fj7j6kzow4vpbxE\nlssKXtASdLryFylEtm/3Jo57/30YOxa6dvW6jWa3+LvFxM+PJ8tlkdgukVZntwpOsHLS1NtHRP5m\n2TJvkFhUlDdKuGnTv+/PclnMXTuX4QuGE1MhhvFtx1Pv9HpHfzMJWWr2EZG/adYMli6FO+7w1hHo\n2RN++eWv/REWwQ3n3UBqn1Ta1mhLq+da0WNeD37c82PwgpaAUvIXKaQiIrypIVJToWRJb7roRx/1\n5gz6Q9Gootzd4m429NtA2eJlqf94fUYtGsWeQ3uCF7gEhJp9RMLE2rVer6BffvFOApdd9s8ym3/f\nzMhFI5n/zXxGXTqKO5rcQXRkdOCDFb+ozV9E/OIcvPYa3HMPXHghJCbCWWf9s1zyT8nEz49ny+4t\njG87nitrX6npIkKQkr+I5Mr+/TBuHMyY4U0Xcc89UKzY38s45/jg6w+Inx9PmWJlSGyXSPMqzY/+\nhhIUSv4ickK++cZL/Ckp3qyhnTp500hkl5mVyXOrnmPUolFcdNZFjG0zlnPKnhOcgOVvlPxF5KR8\n+KF3P6BGDZgyBWrV+meZ/en7eeTzR3hk6SPcXP9mRl42kvIlygc+WPmTunqKyElp3x5Wr4bWreGi\ni2DYMNhzRIefEtEluPfSe1nXZx2ZLpO6M+oy/tPxHEg/EJyg5aT4lfzNLM7M0sxsg5kNPcr+08xs\nnpmtNLMUM+uWbd8A32spZtY/D2MXkTxUpIg3NURKCvz4I9StC7Nn/zVh3B9OL3k60y+fzpLblrDs\nx2XUnl6b51c9r+kiCpgcm33MLALYALQBfgSWAzc659KylRkOnOacG25m5YH1wBlAbeBl4AIgA3gf\nuNM5981RPkfNPiIh5LPPvFHCJUt6XUMbNTp6uSXfLyF+fjz70/eT2C6Rdue0C2ygYSy/m32aARud\nc5udc+nAHOCqI8o44FTf81OBnc65DKAu8IVz7pBzLhP4GLj2RAIVkcC66CJYvtybH6hDB+jdG3bu\n/Ge5llVbsuS2JYy8dCS93+tN3ItxrP55deADllzxJ/lXBn7Itr3F91p204EYM/sRWAUM8L2+BrjE\nzMqYWQngcuAovYpFJBRFRnpTQ6SmeiOGY2LgiScgM/Pv5cyM62KuY13vdXSq1Yn2L7Sn+1vd2bJ7\nS3AClxxF5dH7dACSnXOtzewcYL6ZNXDOpZnZeGA+sBdIBjKP9SYJCQl/Po+NjSU2NjaPwhORk1G2\nLEyfDj16eL2CnnzSawq6+OK/l4uOjKZvs750bdCVCUsm0PCJhvRq2ouhLYdSqlip4ARfiCQlJZGU\nlJQn7+VPm38LIME5F+fbHgY459z4bGXeAcY655b4thcAQ51zXx7xXg8BPzjnnjjK56jNX6QAcA5e\neQXi470pIiZMgDPPPHrZLbu3MHLRSN7b+B73XXIfvc7vRZHIIoENuBDL7zb/5UBNM6tmZkWAG4F5\nR5TZDLT1BXMGUAv4xrddwfdvVeAa4KUTCVREQoMZ3Hij1xRUtaq3nvCECXD48D/LVjmtCs9e9Szz\nu87n3Y3vUu+xery67lV0oRd8fg3yMrM4YCreyeIZ59w4M+uF9wtgpplVAmYBlXxVxjrnXvbV/Rgo\nC6QDA51zScf4DF35ixRAmzbB3XfDxo3eALGOHY9ddv7X8xny0RCKRxUnsV0iLau2DFyghZBG+IpI\n0L37rncSqFvXmyrinGPMAJHlspi9ejb3LbqPppWaMq7tOGqVO8qQYsmRRviKSNBdcQWsWeN1EW3e\nHO67D/bt+2e5CIuga8OupPVJo0WVFrT8b0v6vNuH7fu2Bz7oMKbkLyJ5pmhRb2qIlSu9SePq1vVu\nDh/tR33x6OIMaTmE1D6pREdGEzMjhoc+foj96fsDH3gYUrOPiOSbjz+Gfv28rqLTpkH9+scu+/Wv\nXzNi4QiWfL+EB1o9wK0NbyUyIjJwwRZAavMXkZCVkQEzZ0JCgtdLaPRoKFPm2OW/2PIFg+cP5veD\nvzOh7QTiasZpIZljUPIXkZC3Y4d3H+DNN2HMGOje3RtBfDTOOeatn8fQj4ZS+bTKJLZLpEmlJoEN\nuABQ8heRAmPFCq8p6NAhb9RwixbHLpuemc4zyc8wevFo2tZoy5hWY6hWulrggg1x6u0jIgVGkybw\n6acwYABcdx106wbbth29bHRkNHeefycb+m6gRukaNJnZhCHzh/Dbgd8CGnNhpOQvIgFn5s0WmpYG\np58O553ayerMAAAOLUlEQVQHkydDevrRy59a9FRGtxpNyl0p/H7wd2pPr80jnz/CoYxDgQ28EFGz\nj4gE3fr13i+B77+HqVOhXQ5LAqzdvpZhC4axdvtaHm7zMDfUuyEsbwqrzV9ECjznYN48GDjQWzhm\n8mSoXv34dRZ9u4j4+fFEWAQT20/k0mqXBiTWUKE2fxEp8Mzgqqtg3TrvvkDTpl730APHWSK41dmt\nWNZjGQNbDOTWN2/lypevJPWX1IDFXJAp+YtISClWzOsSmpzsnQjq1oXXXz/6KGHwpovoUr8LaX3S\nuKzaZVw661J6vd2LbXuPcRdZADX7iEiIW7jQW0CmUiXvfkBMzPHL/3rgVx7+5GGeXfks/Zr1Y/BF\ngzmlyCmBCTbA1OwjIoVW69ber4BOnbzFYwYNgl27jl2+bPGyTGw/kS97fMnGXzdS69FazPxqJhlZ\nGYELugBQ8heRkBcd7fUGWrsWdu/2moJmzYKsrGPXObvM2cy+djbzuszj5TUv0+DxBry9/m0tJOOj\nZh8RKXCWL4e+fb2bxI8+ChdccPzyzjne2/geQz4aQoUSFUhsl8gFlXOoVACoq6eIhJ2sLHjuORgx\nwltL4OGHvQFjx5ORlcGslbO4P+l+Lq12KQ+1fogaZWoEJuB8oDZ/EQk7ERHe5HBpaXDqqVCvnjdt\ndMZxmvajIqK4o8kdbOi7gZjyMTR7qhmDPhjErwd+DVzgIULJX0QKtFKlvGUjFy/2Bok1bgyLFh2/\nTskiJRl52UjW9l7LwYyD1J5em8QliRzMOBiYoEOAmn1EpNBwzhsTcM890KwZTJwIVavmXC9tRxrD\nPhpG8rZkHmr9EDfVv4kIC/1rY7X5i4hks38/jB/vTRk9cCAMHuwNHsvJJ5s/IX5+PIczD5PYLpE2\nNdrkf7AnQclfROQovv3WGxewerXXNPSvf3k9hI7HOcfcdXMZvmA4tcvVZkK7CZx3+nmBCTiXlPxF\nRI7jww+9cQLVq8OUKVC7ds51Dmce5okvn+ChTx6i07mdeKDVA1Q+rXK+x5ob6u0jInIc7dvDqlXQ\nti20bAlDhsCePcevUySyCP2b92d93/VUKFmBBk804L6F97H70O7ABJ3PlPxFJCwUKeLdCE5Jge3b\noU4dePHFY08Y94fSxUozru04knsl88PuH6g9vTaPLX+M9MxjrDxTQKjZR0TC0uefe2sJFyvmjRJu\n3Ni/eiu3rWTI/CFs3rWZcW3GcXWdq4O2kIza/EVETkBmJvz3vzByJFxzDYwZA+XK+Vf3g00fED8/\nntOKnkZiu0QuPOvC/A32KNTmLyJyAiIjoUcPSE2FqChvwrjHHvNOCjnpULMDyb2SuaPJHfz71X/T\n+X+d2bhzY/4HnUf8Sv5mFmdmaWa2wcyGHmX/aWY2z8xWmlmKmXXLtm+gma0xs9VmNtvMiuRh/CIi\nJ61MGa/p56OP4H//81YR++STnOtFRkTSrVE31vddT9NKTbnwmQvp/35/ftn3S/4HfZJyTP5mFgFM\nBzoA9YAuZlbniGJ9gLXOuUZAK2CSmUWZ2ZlAP6CJc64BEAXcmJdfQEQkrzRo4E0NMXw43HST99i6\nNed6JaJLMPyS4aT28ZaQrDujLmM/GcuB9OOsQRlk/lz5NwM2Ouc2O+fSgTnAVUeUccCpvuenAjud\nc39MrxQJlDSzKKAE8OPJhy0ikj/M4IYbvAnjzj4bGjaEcePg0KGc61YoWYFpHafx+e2fs2LbCmpN\nr8WslbPIzPKjHSnA/En+lYEfsm1v8b2W3XQgxsx+BFYBAwCccz8Ck4Dvga3A7865j042aBGR/Fay\nJDz0ECxdCkuWQP368N57/tU9t9y5zL1+Lv/r/D+eWvEUTWY24YNNH+RvwLmUVzd8OwDJzrkzgcbA\nDDM7xcxK4/1KqAacCZxiZjfl0WeKiOS7mjXh7be9kcEDBnhTRGza5F/dC8+6kE+7f0rCZQn0e78f\n7V9oz6ptq/I3YD9F+VFmK5B9Xrwqvtey6w6MBXDOfW1m3wJ1gOrAN865XwHM7HXgIuClo31QQkLC\nn89jY2OJjY31IzwRkfx3+eXQpo13EmjRAnr2hHvv9X4hHI+ZcU3da+hUqxNPrXiKDi92oEPNDoxp\nNYazSp2VqxiSkpJISko68S+RPa6c+tabWSSwHmgD/AQsA7o451KzlZkBbHfOjTazM4AvgYZATeAZ\n4ALgEPAssNw5N+Mon6N+/iJSIGzd6k0R8fHHkJjo3SPwd5zX7kO7mbBkAo9/+Tg9mvRg+MXDKVWs\n1AnFke+DvMwsDpiK10z0jHNunJn1ApxzbqaZVQJmAZV8VcY651721b0fr4dPOpAM3OG7cXzkZyj5\ni0iB8skn3ijh0qW9VcQaNPC/7tbdWxm1aBTvbHyHEReP4K4L7qJIZO56wmuEr4hIkGRmwpNPQkKC\n9wtg9GgoW9b/+ik/pzD0o6Fs2LmBsW3G0jmms9/TRSj5i4gE2c6dcN993kpiDz4It9/ujSD214Jv\nFhA/P54ikUWY2H4iF1e9OMc6Sv4iIiEiOdlrCjp40Bs1fGEupvzJclm8lPIS9y68l8YVGzO+7Xhq\nlz/24gOa20dEJEQ0buzdC7j7bujcGW69FbZt869uhEVwS4NbWN93PS3PasnFz15M73d78/Pen/M8\nTiV/EZE8Zga33OKNEq5YEc47z1tM/vBh/+oXiypGfMt40vqkUSyqGPUeq8eDix9k3+F9eRajkr+I\nSD459VRvIfklS2DBAq830Icf+l+/XIlyTO4wmWU9lrFuxzpqTa/F0yueJiMrI+fKOVCbv4hIADjn\njRQeONA7CUye7M0dlBvLti4jfn48O/fvZHzb8XSq3Uk3fEVECoKDB2HSJC/59+kDw4ZBiRL+13fO\n8c6Gdxjy0RDS+qYp+YuIFCQ//ACDB3sTx02aBNdd5/8oYYCMrAyiI6OV/EVECqJFi6B/fzj9dG+U\ncL16/tdVV08RkQKqVStvbMBVV0FsrHdP4Pff8/9zlfxFRIIsKsq7+l+7Fvbu9dYS/u9/ISsr/z5T\nzT4iIiFm+XJvlLBz3ijhZs2OXk7NPiIihcgFF8Bnn8Fdd3nNQbffDtu35+1nKPmLiISgiAjo1s0b\nJVyqlHcjeOpUSP/HhPgnRs0+IiIFwLp13jKSP/3k9Qpq3VqzeoqIhAXn4I03YNAgr2no1VfV5i8i\nUuiZwbXXer8CcjMe4KjvFSpX27ryFxHJHfX2ERGRXFHyFxEJQ0r+IiJhSMlfRCQMKfmLiIQhJX8R\nkTCk5C8iEoaU/EVEwpCSv4hIGFLyFxEJQ34lfzOLM7M0M9tgZkOPsv80M5tnZivNLMXMuvler2Vm\nyWa2wvfvLjPrn8ffQUREcinH5G9mEcB0oANQD+hiZnWOKNYHWOucawS0AiaZWZRzboNzrrFzrgnQ\nFNgHvJGn36AQSkpKCnYIIUHH4S86Fn/Rscgb/lz5NwM2Ouc2O+fSgTnAVUeUccCpvuenAjudcxlH\nlGkLfO2c++FkAg4H+uP26Dj8RcfiLzoWecOf5F8ZyJ6wt/hey246EGNmPwKrgAFHeZ8bgJdPJEgR\nEclbeXXDtwOQ7Jw7E2gMzDCzU/7YaWbRwJXA3Dz6PBERORnOueM+gBbA/2XbHgYMPaLMO0DLbNsL\ngPOzbV+Z/T2O8TlODz300EOP3D1yyuHHekSRs+VATTOrBvwE3Ah0OaLMZrw2/SVmdgZQC/gm2/4u\n5NDkc6ILEoiISO75tZKXmcUBU/GaiZ5xzo0zs154Z52ZZlYJmAVU8lUZ65x72Ve3BN7JoYZzbk8+\nfAcREcmlkFnGUUREAiegI3xzGizmKzPNzDb6Bow1CmR8geTHwLmbzGyV7/GpmdUPRpyB4M/fha/c\nBWaWbmbXBjK+QPLz/0isb9DkGjNbFOgYA+VEB5cWRmb2jJn9bGarj1Mmd7nzRG8W5PaBd6LZBFQD\nooGVQJ0jynQE3vU9bw4sDVR8gXz4eSxaAKV8z+PC+VhkK7cAr3PBtcGOO4h/F6WAtUBl33b5YMcd\nxGMxHK+JGaA8sBOICnbs+XQ8LgYaAauPsT/XuTOQV/7+DBa7CngewDn3BVDKdwO5sMnxWDjnljrn\ndvk2l/LPsRWFhT9/FwD9gFeB7YEMLsD8ORY3Aa8557YCOOd2BDjGQMmrwaWFgnPuU+C34xTJde4M\nZPL3Z7DYkWW2HqVMYeDPscjuDuD9fI0oeHI8FmZ2JnC1c+5xoDD3CvPn76IWUNbMFpnZcjPrGrDo\nAiuvBpeGi1znTn+6ekoQmVkroDvez75wNQXI3uZbmE8AOYkCmgCtgZLA52b2uXNuU3DDCoo/Bpe2\nNrNzgPlm1sA5tzfYgRUEgUz+W4Gq2bar+F47ssxZOZQpDPw5FphZA2AmEOecO95PvoLMn2NxPjDH\nzAyvbbejmaU75+YFKMZA8edYbAF2OOcOAgfN7GOgIV77eGHiz7HoDowFcM59bWbfAnWALwMSYWjJ\nde4MZLPPn4PFzKwI3mCxI//zzgP+A2BmLYDfnXM/BzDGQMnxWJhZVeA1oKtz7usgxBgoOR4L51wN\n3+NsvHb/3oUw8YN//0feAi42s0jfGJrmQGqA4wwEf47FH4NLOcbg0sLGOPav3lznzoBd+TvnMs2s\nL/Ahfw0WS80+WMw5956ZXW5mm/Cmf+4eqPgCyZ9jAYwEygKP+a54051zzYIXdf7w81j8rUrAgwwQ\nP/+PpJnZB8BqIBOY6ZxbF8Sw84WffxdjgFnZuj8Occ79GqSQ85WZvQTEAuXM7HvgfqAIJ5E7NchL\nRCQMaRlHEZEwpOQvIhKGlPxFRMKQkr+ISBhS8hcRCUNK/iIiYUjJX0QkDCn5i4iEof8H+sQDGrO7\nCCEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x90c88d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alp = [0.001, 0.003, 0.01, 0.03, 0.1, 0.3, 1]\n",
    "plt.plot(alp, accuracy_tra, color=\"blue\", linewidth=1.0, linestyle=\"-\")\n",
    "plt.plot(alp, accuracy_val, color=\"green\", linewidth=1.0, linestyle=\"-\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们可以看到学习率为0.01时，在验证集上的表现最好，让我们以该学习率在训练集上再训练一次："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "------------- Training model for number: 0 -------------\n",
      "\n",
      "step 0: 0.049542\n",
      "step 10: 0.030585\n",
      "step 20: 0.027832\n",
      "step 30: 0.026479\n",
      "step 40: 0.025645\n",
      "step 50: 0.025066\n",
      "step 60: 0.024633\n",
      "step 70: 0.024291\n",
      "step 80: 0.024010\n",
      "step 90: 0.023773\n",
      "Final training error: 0.023587\n",
      "\n",
      "------------- Training model for number: 1 -------------\n",
      "\n",
      "step 0: 0.046138\n",
      "step 10: 0.032117\n",
      "step 20: 0.029923\n",
      "step 30: 0.028764\n",
      "step 40: 0.028000\n",
      "step 50: 0.027442\n",
      "step 60: 0.027007\n",
      "step 70: 0.026653\n",
      "step 80: 0.026356\n",
      "step 90: 0.026100\n",
      "Final training error: 0.025898\n",
      "\n",
      "------------- Training model for number: 2 -------------\n",
      "\n",
      "step 0: 0.096797\n",
      "step 10: 0.073934\n",
      "step 20: 0.071076\n",
      "step 30: 0.069663\n",
      "step 40: 0.068759\n",
      "step 50: 0.068106\n",
      "step 60: 0.067601\n",
      "step 70: 0.067192\n",
      "step 80: 0.066851\n",
      "step 90: 0.066561\n",
      "Final training error: 0.066334\n",
      "\n",
      "------------- Training model for number: 3 -------------\n",
      "\n",
      "step 0: 0.112371\n",
      "step 10: 0.087489\n",
      "step 20: 0.084005\n",
      "step 30: 0.082545\n",
      "step 40: 0.081721\n",
      "step 50: 0.081171\n",
      "step 60: 0.080766\n",
      "step 70: 0.080447\n",
      "step 80: 0.080187\n",
      "step 90: 0.079968\n",
      "Final training error: 0.079798\n",
      "\n",
      "------------- Training model for number: 4 -------------\n",
      "\n",
      "step 0: 0.083006\n",
      "step 10: 0.058826\n",
      "step 20: 0.055567\n",
      "step 30: 0.053897\n",
      "step 40: 0.052806\n",
      "step 50: 0.052002\n",
      "step 60: 0.051371\n",
      "step 70: 0.050854\n",
      "step 80: 0.050420\n",
      "step 90: 0.050048\n",
      "Final training error: 0.049755\n",
      "\n",
      "------------- Training model for number: 5 -------------\n",
      "\n",
      "step 0: 0.125051\n",
      "step 10: 0.094596\n",
      "step 20: 0.089975\n",
      "step 30: 0.087658\n",
      "step 40: 0.086194\n",
      "step 50: 0.085156\n",
      "step 60: 0.084366\n",
      "step 70: 0.083736\n",
      "step 80: 0.083215\n",
      "step 90: 0.082774\n",
      "Final training error: 0.082430\n",
      "\n",
      "------------- Training model for number: 6 -------------\n",
      "\n",
      "step 0: 0.062390\n",
      "step 10: 0.043676\n",
      "step 20: 0.041088\n",
      "step 30: 0.039839\n",
      "step 40: 0.039078\n",
      "step 50: 0.038553\n",
      "step 60: 0.038161\n",
      "step 70: 0.037851\n",
      "step 80: 0.037597\n",
      "step 90: 0.037381\n",
      "Final training error: 0.037213\n",
      "\n",
      "------------- Training model for number: 7 -------------\n",
      "\n",
      "step 0: 0.068695\n",
      "step 10: 0.052122\n",
      "step 20: 0.050118\n",
      "step 30: 0.049062\n",
      "step 40: 0.048355\n",
      "step 50: 0.047828\n",
      "step 60: 0.047411\n",
      "step 70: 0.047068\n",
      "step 80: 0.046778\n",
      "step 90: 0.046528\n",
      "Final training error: 0.046329\n",
      "\n",
      "------------- Training model for number: 8 -------------\n",
      "\n",
      "step 0: 0.168364\n",
      "step 10: 0.121003\n",
      "step 20: 0.116161\n",
      "step 30: 0.114525\n",
      "step 40: 0.113700\n",
      "step 50: 0.113178\n",
      "step 60: 0.112803\n",
      "step 70: 0.112514\n",
      "step 80: 0.112280\n",
      "step 90: 0.112086\n",
      "Final training error: 0.111936\n",
      "\n",
      "------------- Training model for number: 9 -------------\n",
      "\n",
      "step 0: 0.136421\n",
      "step 10: 0.105690\n",
      "step 20: 0.101991\n",
      "step 30: 0.100363\n",
      "step 40: 0.099392\n",
      "step 50: 0.098722\n",
      "step 60: 0.098221\n",
      "step 70: 0.097827\n",
      "step 80: 0.097507\n",
      "step 90: 0.097239\n",
      "Final training error: 0.097032\n",
      "\n",
      "training time 594.169s\n"
     ]
    }
   ],
   "source": [
    "from time import time\n",
    "\n",
    "t0 = time()\n",
    "predictors = list()\n",
    "for number in range(0, PREDICT_NUMBERS):\n",
    "    train_y = np.array(list(map(int, pd_train[0].values == number)))\n",
    "    predictors.append(LogisticRegression())\n",
    "    print(\"------------- Training model for number: %d -------------\\n\" % number)\n",
    "    predictors[number].fit(train_X, train_y, 0.01, batch_num=10, verbose=True)\n",
    "    \n",
    "print('training time %.3fs' % round(time() - t0, 3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "训练时间564s\n",
    "\n",
    "### 5.测试\n",
    "\n",
    "最后我们把训练好的模型，用在测试集上，看看准确率如何："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Test model on testing data\n",
      "\n",
      "On testting data: 9183/10000\n",
      "\n",
      "testting time: 2.476s\n"
     ]
    }
   ],
   "source": [
    "t0 = time()\n",
    "\n",
    "print(\"Test model on testing data\")\n",
    "test_X = pd_test.drop([0], axis=1).values / 255.\n",
    "tot_right_num = 0\n",
    "for tx, ty in zip(test_X, pd_test[0].values):\n",
    "    pred_y = list()\n",
    "    for number in range(0, PREDICT_NUMBERS):\n",
    "        pred_y.append(predictors[number].predict(np.reshape(tx, (1, -1))))\n",
    "\n",
    "    if np.argmax(pred_y) == ty:\n",
    "        tot_right_num += 1\n",
    "\n",
    "print(\"\\nOn testting data: %s/%s\\n\" % (tot_right_num, test_X.shape[0]))\n",
    "\n",
    "print(\"testting time: %.3fs\" % round(time() - t0, 3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "最终我们在测试集上的准确率为91.8%，比起baseline有所提高；在10000个样本上的测试时间为2.48s，这个结果是令人满意。\n",
    "\n",
    "### 6.持久化\n",
    "\n",
    "我们把训练好的模型储存起来："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import pickle\n",
    "\n",
    "with open('logistic_regression.pkl', 'wb') as f:\n",
    "    pickle.dump(predictors, f)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "测试持久化是否成功："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Test model on testing data\n",
      "\n",
      "On testting data: 9183/10000\n",
      "\n"
     ]
    }
   ],
   "source": [
    "clf = pickle.load(open('logistic_regression.pkl', 'rb'))\n",
    "\n",
    "print(\"Test model on testing data\")\n",
    "test_X = pd_test.drop([0], axis=1).values / 255.\n",
    "tot_right_num = 0\n",
    "for tx, ty in zip(test_X, pd_test[0].values):\n",
    "    pred_y = list()\n",
    "    for number in range(0, PREDICT_NUMBERS):\n",
    "        pred_y.append(clf[number].predict(np.reshape(tx, (1, -1))))\n",
    "\n",
    "    if np.argmax(pred_y) == ty:\n",
    "        tot_right_num += 1\n",
    "\n",
    "print(\"\\nOn testting data: %s/%s\\n\" % (tot_right_num, test_X.shape[0]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们可以发现，从刚才生成的**logistic_regression.pkl**文件中读取的clf一样可以进行预测。\n",
    "\n",
    "### 7.总结\n",
    "\n",
    "在项目中，我们用Logistic regression算法实现了一个手写识别系统，最终在测试集上的准确率为91.8%。\n",
    "\n",
    "学习了：\n",
    "\n",
    "-  Logistic regression用于多分类问题\n",
    "- 学习速率的选择\n",
    "- 使用了pickle进行持久化"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python [Root]",
   "language": "python",
   "name": "Python [Root]"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
